Semiclassical Gaps in the Density of States of Chaotic Andreev Billiards 
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The connection of a superconductor to a chaotic ballistic quantum dot leads to interesting phe- 
nomena, most notably the appearance of a hard gap in its excitation spectrum. Here we treat such 
an Andreev billiard semiclassically where the density of states is expressed in terms of the classical 
trajectories of electrons (and holes) that leave and return to the superconductor. We show how 
classical orbit correlations lead to the formation of the hard gap, as predicted by random matrix 
theory in the limit of negligible Ehrenfest time Te, and how the influence of a finite te causes the 
gap to shrink. Furthermore, for intermediate te we predict a second gap below E = 7t/i/2te which 
would presumably be the clearest signature yet of TE-effects. 
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A superconductor (S) in contact with a normal con- 
ductor (N) considerably affects its spectral density of 
quasiparticle excitations: due to Andreev reflection [l[ 
at the NS interface the density of states (DoS) is sup- 
pressed closely above the Fermi energy E-p. This prox- 
imity effect is also expected for an 'Andreev billiard' 0, 
an impurity-free quantum dot attached to a supercon- 
ductor 0, |j| , and has attracted considerable theoretical 
attention during the last decade (see @ for a review). 

An Andreev billiard has the interesting peculiarity that 
the suppression of its (mean) DoS crucially depends on 
whether the dynamics of its classical counterpart is in- 
tegrable or chaotic: while the DoS vanishes linearly in 
energy for the intcgrable case, the spectrum of a chaotic 
billiard is expected to exhibit a true gap above E-p 
Based on random matrix theory (RMT) this gap was 
predicted to scale with the Thouless energy, Et = H/2te>, 
where td is the average (classical) dwell time a particle 
stays in the billiard between successive Andreev reflec- 
tions 0]. On the contrary, semiclassics based on the so- 
called Bohr-Sommerfeld (BS) approximation yields only 
an exponential suppression of the DoS @-@], a discrep- 
ancy that has attracted much theoretical interest jl0l - ll4| . 
Lodder and Nazarov Q pointed out that these seem- 
ingly contradictory predictions are valid in different lim- 
its, governed by the ratio r = te /td . Here the (quantum 
mechanical) Ehrenfest time te ~ | In H\ separates the evo- 
lution of wave packets following essentially the classical 
dynamics from longer time scales dominated by wave in- 
terference. In the universal regime, t = 0, the Thouless 
gap (from RMT) is clearly established [f| , while the 
BS approximation describes the classical limit t— >oo. 

Various approaches have been used to better under- 
stand the crossover from the Thouless to the Ehrenfest 
regime of large t, where RMT loses its applicability fioj ]. 
These include effective RMT [12(, predicting a gap size 
scaling with the Ehrenfest energy Ep = h/2rp, as well as 
stochastic [HI and perturbative 11 1 methods. Recently 
the gap at irEp was confirmed for t> 1 in a quasiclassical 



approach based on the Eilenberger equation [f3 |. 

The purpose of this Letter is twofold. Firstly, using the 
scattering approach (IBj . we demonstrate that the DoS 
can be evaluated semiclassically for te = by using an 
energy-dependent extension of the work flil ] on the mo- 
ments of the transmission eigenvalues. This semiclassi- 
cally computed DoS yields a hard gap, in agreement with 
RMT. Secondly we address the whole crossover regime of 
t > 0, by incorporating the te dependence. In the limit 
r>l, the width of the gap approaches 7t.Ee, eventually 
recovering the BS prediction for t— > oo. More interest- 
ingly in the intermediate regime t > 1 we predict the 
appearance of a second 'Ehrenfest' gap at irEp. 

Andreev billiard. — In the scattering approach the su- 
perconductor is represented by a lead that carries M 
scattering channels, and the excitation spectrum can be 
entirely expressed in terms of the (electron) scattering 
matrix S [HI]. The average DoS reads (when divided 
by twice the average density of the isolated billiard), 

„ (e) = 1+2 £fc£ In W, (1) 

n ae 

n—l 

in terms of correlation functions of n S'-matrices, 



C(e,n) = — Tr 



, (2) 



at different energies. Here the energy difference 2E is ex- 
pressed in units of the Thouless energy and td = Tn/M 
with Th the Hcisenbcrg time, i.e. the time conjugate to 
the mean level spacing. For e = 0, the C(e, n) in Eq. ((2J, 
with transmission rather than scattering matrices, give 
the moments of the transmission eigenvalues which were 
calculated semiclassically (to leading order in Af _1 ) in 
[lit and our derivation is based on that work. 

Semiclassical evaluation in the universal regime. — To 
evaluate Eq. ([2]) , we start with the semiclassical approxi- 
mation to the scattering matrix elements connecting the 
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channel a to b, which are given by [Tj 



Sba(EF ± E) 
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in terms of the classical trajectories £ connecting a to fe. 
Here is the action of and is its stability (includ- 
ing Maslov indices). We substitute Eq. §5§ into Eq. J2]) 
and expand the action up to first order in the energy 
yielding the duration Tq = dS^/dE-p. The correlators 
are then given by a sum over 2n trajectories 

1 " 

c(e ' n)sa M7^nE E E 

j = 1 a-j ,bj Cj ( aj -> b j ) q ( b 3 ; -> a 3 + j ) 

x exp [(i/S)(5 c . -S c ,)] exp [(ie/2r D )(T Cj + T c; )], (4) 

with a„-|-i = oi. The final trace of the product of matri- 
ces means that the trajectories complete a cycle, moving 
forward along the Q and back along the ; for an exam- 
ple of this structure for n = 3, see Fig. [TJa). 

In Eq. ((4]) we add the actions of all the unprimed tra- 
jectories and subtract the actions of the primed ones, 
so the resulting phase oscillates wildly, unless the total 
action difference is of the order of h. One way to get 
small action differences is to collapse all the trajectories 
onto each other [see Fig. [T{b)]. This leads to encounters 
where the electron trajectories avoid crossing while the 
hole trajectories cross (or vice versa) to ensure that they 
each connect the correct channels. In phase space, the 
trajectories of course do not cross, but just come close 
enough together to allow this reconnection. Besides this 
direct collapse further possibilities arise from sliding the 
encounters together or into the leads (see (lj|). 

For each possibility we also need to know its semi- 
classical contribution. Following the treatment for open 
systems of the first off- diag onal pair by [18|, the gen- 
eralization to all orders [19( led to diagrammatic rules, 
whereby each link (i.e. each trajectory stretch connect- 
ing channels or encounters) essentially gives a factor of 
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Figure 1: (a) Schematic picture of the trajectory structures 
for n — 3. The (blue) solid lines represent electrons which are 
retro-reflected as holes [dashed (red) lines], (b) Collapsing 
the trajectories onto each other leads to encounters, (c) The 
end result, i.e. correlated Andreev reflected paths. 



[M(l — ie)] -1 , while each i-encounter (where I electron 
trajectories come close together) contributes —M(l —He) 
as the encounter stretches all remain inside the cavity or 
touch the lead together. Summing the contributions, by 
extending the work of [lj| to include energy differences 
[20| . and using the diagrammatic rules above, we arrive 
at the intermediate generating function g(e, r), which in- 
cludes all possible diagrams apart from where the top 
encounter enters the lead, and is given implicitly by 



g (1 - ie) - 1 = rg 2 (g - 1 - ie) 



(5) 



Including the possibility where the top encounter can en- 
ter the lead, leads to the generating function 



G(e,r) 



71=1 
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1 - rg 



(6) 



of the correlation functions. By inverting Eq. Q we can 
see that G is given implicitly by the cubic equation 



r(r-iyG 6 + r{3r + ie-3)G z + (3r + ie-l)G = -1. (7) 

Expanding G (or g) as a power series in r, we obtain the 
first couple of correlation functions (which can be checked 
by considering the semiclassical diagrams explicitly) as: 



C(e,l) 



II -ie)' 



CM) = iz", (8) 

(1 - ie) 



Density of states in the universal regime. — We can gen- 
erate these correlation functions recursively to obtain a 
truncated version of the sum in Eq. (fl}. However, we can 
go one step further and find the generating function of 
the terms that appear in the density of states 



H(e,r) 



fr de 



G(e,r) dr 



~ r n-l QC^ n} 



71=1 



dc 



d(e) = l-2Re#(e,-l), 



(9) 



using a computer aided search over cubic equations with 
low order polynomial coefficients to obtain 

(er) 2 (l - r)H 3 + ier [r(ie - 2) + 2(1 - ie)] H 2 

+ [r(l - 2ie) - (1 - ie) 2 ] H + 1 = 0. (10) 

As G and H are solutions of algebraic generating func- 
tions, so too must be their derivatives, and we prove the 
result in Eq. (fTOj) by differentiating G with respect to e 
and rH with respect to r and demonstrating that these 
derivatives satisfy the same algebraic equation. 

Taking the solution whose expansion agrees with 
Eq. (g]), the DoS then follows from Eq. © as 

d(e) = V3/(6e) [Q+(e) - Q_(e)] , e > 2[(V5 - l)/2]« , 

(11) 

where Q±(e) = [8 - 36e 2 ± 3eV3e 4 + 132e 2 - 48J ' . This 
is exactly the RMT result @] (dotted line in Fig.©. 
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Density of states in the Ehrenfest regime. — The effect 
of non-zero Ehrenfest time te on the first three correla- 
tion functions (7(e, r, n) with r = te/td has previously 
been calculated semiclassically [2l[ . For these, the effect 
of increasing te is twofold; first as each encounter typ- 
ically lasts te, forming the diagrams considered before 
becomes less likely, while conversely the possibility that 
all the trajectories are correlated for their whole length 
increases (cf the bands in |l_4|), and we need to add this 
additional set of diagrams. For n<3 the results 2l| sug- 
gest the replacement (which is in line with effective RMT 



C(e,T,n) = C(e,n)c- r(1 ~ in£) + 



1 _ g-T(l-ine) 

1 — me 



(12) 



This separation into two terms was shown in [22| (for 
e = 0) to be a consequence of the preservation under time 
evolution of phase-space volume inside the system and 
hence the splitting of phase-space into two essentially in- 
dependent subsystems. Including an energy difference, 
the second term in Eq. ([T2|) as well as the form of the 
exponential in the first term follow directly. The only 
thing that cannot be determined from [22| is the remain- 
ing factor in the first term. To show that it is indeed 
C(e, n) we reconsider the diagrams treated before, which 
were created by sliding the encounters together or into 
the lead, as part of a continuous deformation of a single 
family of diagrams. By suitably transforming their semi- 
classical contributions, we can extract their complete te 
dependence and show it is always e~ T ( 1 ~ me \ Summing 
over all the families, this common te dependence, along 
with the necessity of recovering the previous result for 
t = 0, shows that Eq. (fT2j) holds for all n [H]. It is possi- 
ble to obtain the rest of this equation semiclassically for 
all n. By considering correlated trajectory bands explic- 
itly we obtain the second term, while the separation into 
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Figure 2: Density of states for t = te/td ~2 (solid line), along 
with the BS (dashed) and RMT (dotted) limits, showing a 
second gap just below £ = ■§■ Inset: TE-related 7r-periodic DoS 
oscillations at higher energy (after subtracting the BS curve). 



two terms arises from opposing restrictions on the bands 
and the encounters p3| . 

Equation (fT2"j) reproduces the two known limits: the 
previous RMT result for t = 0, and the BS result 

<fes(e) = (7) s°n$(J/% ' for T = 00 ' Alon S side the two 
limits, this equation, and in particular the ne depen- 
dence in the exponent, leads to interesting TE-effects: a 
re-normalized gap, an oscillatory DoS with spikes with 
period 27r/r, and a second intermediate gap. 

To study this behavior, we substitute Eq. (|12p into 
Eq. ([T]) and get two contributions, from the two terms. 
The first yields a reduced RMT-type contribution that 
can be evaluated as before. The energy differential leads 
to two further terms, and this part of the DoS is 

dj(e) = e~ T [1 - 2Re e ieT H(e, -e ier )] (13) 
+ rc- r [1 - 2Rc e ier G(e, -e ieT )] , 

which includes the most natural constant term and re- 
duces to © when r = 0. Including the rest of the con- 
stant term (=1 — (1 + r)e~ r ) from Eq. ([1]) with the second 
contribution of (|12[) . which can be summed exactly via 
Poisson summation, we obtain 



00 f—Fi™ fl 
(i 2 ( £ ) = l-(l+r)c- T +2y^Im 

Z ^ 77 



= cfes(e) - e" 



n=l 



efes (e) 



dc 



2_ c -r(l-i«£) 

1 — me 



2k (ir/ef 
sinh (ir/e) 



(14) 



where k = [ £ T ^ \ involving the floor function. We note 
that di{e) is zero up to er = 7r, but when combined with 
the first contribution, we find interesting new features. 

As an illustration we plot the full DoS for t = 2 in 
Fig. [2] We find a clear reduction of the RMT gap and 
in the inset an oscillatory behavior of the DoS at larger 
energy. We note that TE-oscillations have previously been 
predicted 11, HI, however those appearing here have a 
larger magnitude. More interestingly though, the result 
in Fig. [5] shows the appearance of a second pronounced 
gap. This structure in the DoS would be a clear-cut 
signature of the Ehrenfest time. The absence of such a 
feature in previous numerical work is presumably due to 
the difficulty in reaching the limit r > 1. 

More generally, using Eqs. Q and (|10|) for G and H, 
we can express d\(e) in (|13p explicitly in a form similar 
to Eq. (fTTj) . The result is only non-zero when 



D(e,r) = e 4 -8e 3 sin(er)+4e 2 [5 + 6cos(er)] 
+24esin(er) - 8 [1 + cos(er)] , 



(15) 



is positive. D{e 1 r) is negative up to the first root of 
Eq. (fj"5|) . and so we see a hard gap up to this point. As 
t is increased this first gap shrinks and eventually ap- 
proaches er = 7r, i.e. E = ttEe for t>1. We recall that 
the second contribution to the DoS from Eq. (QJJ is ex- 
actly zero up to this point, so if we take the limit r— >oo 
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Figure 3: (a), (b) Density of states as a function of €t = E/Ee 
for various values of r showing the appearance of a second gap 
below er — n. Inset: DoS for r = 20 (solid line) together with 
the BS limit (dashed), (c) Width of the original (first) gap 
as a function of r. Our semiclassical result (solid line) agrees 
with effective RMT 0] (dotted line); dashed line: prediction 
from the stochastic model of [l3l |. (d) Width of the second 
gap as a function of r. 



is so small as to be essentially indistinguishable from BS, 
apart from perhaps the following spike and second gap. 

Conclusions. — Based on a systematic treatment of cor- 
relation functions involving n scattering matrices, we cal- 
culated the DoS of an Andrcev billiard semiclassically, 
and recover a hard gap extending up to 0.6-Et as in RMT 
(at te = 0). Likewise, increasing te we can see how the 
gap closes (approaching E — ttE-e) in agreement with ef- 
fective RMT, and we can study the full crossover from 
the RMT limit to te 3> td ■ Interestingly this transition is 
not smooth, and inbetween we see the formation of a sec- 
ond gap at Ec^ttE-e for te — tq. Such a striking feature, 
which would be interesting to confirm by independent 
means, should be an easier TE-signature to observe ex- 
perimentally than the change in size of the original gap. 

We thank I. Adagideli, A. Altland, Ph. Jacquod, 
M. Novaes, J.D. Urbina, and R.S. Whitney for valuable 
discussions. We acknowledge funding from the DFG un- 
der GRK 638 (DW, KR), the NSF under Grant 0604859 
(GB) and from the AvH Foundation (JK, CP). 



at fixed er, we observe a hard gap up to ttEe in agree- 
ment with the recent complementary quasiclassical work 
of 14|. Away from this limit though, and for general 
t, alongside the zero contribution (up to er = tt) from 
the bands of correlated trajectories, we also have to in- 
clude the contribution from the trajectories with encoun- 
ters which determine the exact size of the first gap and 
cause the behavior described below. We plot the width 
of this gap in Fig. Etc) and find exact agreement with 
the effective RMT prediction [l2[ ; it also seems to be 
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in accordance with previous numerical findings [2 
limited to r < 1 . 

But when r > 0.916, the discriminant (|15|) has addi- 
tional roots. In particular D(e,r) is negative between 
the second and third roots so that the first contribution 
to the DoS again falls to zero and we see the creation of 
a second gap. As r is increased, the roots spread apart 
and the gap widens, but of course G?2(e) is only zero up 
to £T = 7r. For r > ^ the third root of Eq. (fl~5j) is beyond 
this and so the second gap is cut short by the jump at 
er = tt in the DoS coming from the second contribution in 
Eq. (fT4j) . As r is increased further the second gap starts 
to shrink, as can be seen from the plot of its width in 
Fig. Eld). 

In Figs.[3Ja) and[3{b) we show the DoS for different val- 
ues of r, illustrating the formation and then the shrinking 
of the second gap. Though it shrinks, the second gap per- 
sists even for large values of r as can be seen in the inset. 
Also visible in Figs. Eta) and Etb) is that the first gap 
slowly approaches er = tt: for r = 20 the first hard gap 
ends at ex « 2.44, but for such a large r the contribution 
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